25/11/20 PH
Build on VMI class for VMI data analysis:
Implemented:
To do:
Background:
# # Memory profiler - OPTIONAL for testing
# # https://timothymonteath.com/articles/monitoring_memory_usage/
# %load_ext memory_profiler
# %memit
# # Quick hack for slow internet connection!
# %autosave 36000
# Standard imports
from pathlib import Path
import sys
import numpy as np
import xarray as xr
# Import class from file
sys.path.append(r'/cds/home/p/phockett/dev/tmo-dev/tmo')
import tmoDataBase as tb
# from vmi import VMI # VMI class
from inversion import VMIproc # VMI processing class
tb.setPlotDefaults() # Set plot defaults (optional)
# Set main data path
baseDir = Path('/reg/data/ana15/tmo/tmolw0618/scratch/preproc/v2')
# Setup class object and which runs to read.
data = VMIproc(fileBase = baseDir, runList = range(89,97+1), fileSchema='_preproc_elecv2')
# The class has a bunch of dictionaries holding info on the data
# data.runs['files']
# Read in the files with h5py
# There is very basic checking implemented currently, just to see if 'energies' is present.
data.readFiles()
Note this defaults to electron coords ('xc','yc').
data.genVMIXmulti()
# Check an image with Xarray plotter (Matplotlib)
data.imgStack['signal'].sel(run=89).pipe(np.log10).plot.imshow()
This is set as native pixel value, or assumed to be the max point in the image.
data.setCentre([553, 546]) # Set explictly
# data.setCentre() # Default case - set as max ind.
data.checkCentre()
# All the coord parameters are set in `centreInds`
# Note here that `cList` gives the array indicies for the centre point, based on the pixel coords set.
data.centreInds
See the previous demo for info, here we'll just generate subtracted images.
# Subtraction, and set as new image set
data.imgStack['sub'] = data.imgStack['signal'] - data.imgStack['bg']
data.imgStack
# View subtracted image stack
data.showImgSet(name='sub')
Multiple methods to be supported... so far only Elio's cpBasex is implemented, and this needs a bit of setup.
# For cpbasex, set the basis path
basisPath = r'/reg/d/psdm/tmo/tmolw0618/scratch/results/tmo_ana/calc/G_r512_k128_l4.h5'
# Set also local version - pip installed version not working with loadG in testing (version mismatch issue, or something in import routine?)
pbasexPath = r'/reg/d/psdm/tmo/tmolw0618/results/modules/pbasex' #/pbasex'
# Import module + basis functions
data.setMethod(pbasexPath = pbasexPath, basisPath = basisPath)
Invert a set of images by setting filterSet (if not set, all images in self.imgStack will be processed). Optionally, set normalisation parameters for the images, the default is norm={'type':'max', 'scope':'global'}, which normalises to the global maximum (i.e. over all images).
Currently the options are:
type:
scope:
# Invert a dataset
data.inv(filterSet='signal')
# Full results are stored in a dictionay, self.proc[filterSet]
filterSet = 'signal'
data.proc[filterSet].keys()
# ... where 'xr' contains the results in an Xarray
data.proc[filterSet]['xr']
# For rapid plotting per run, use self.plotSpectra()
data.plotSpectra()
# This can also be set to return a HoloMap for more flexibility
hmap = data.plotSpectra(returnMap = True)
hmap.overlay('BLM').layout('run').cols(1)
Various other stages and results are also logged to the data strucure...
Images are set to 'filterSet-TYPE', where TYPE is:
quadrant.unfoldQuadrant for full symmetrized images.)These can be accessed and plotted with the usual functionality.
NOTE: there is currently no default radial masking set, so the default plot scaling will likely be off - adjust with the histogram.
# Quad filter
data.showImg(name='signal-fold')
data.showImg(name='signal-inv')
# Raw cpbasex output data
data.proc['signal']['pbasex'].keys()
A basis Poissonian bootstrapping routine is implemented for estimating uncertainties. This operates on the event data, and allows error propagation through the image generation & analysis, but is (consequently) a little slow. (For more on the method employed, see the wikipedeia page#Poisson_bootstrap).)
Note this routine is a little crude at the moment! It will run image generation for all filterSets, run image subtraction, then image inversion & stats for all filterSets, or just the listed element.
# Run bootstrapping & inversion routine, 5 iterations, and calculate stats for the subtracted data only.
data.bootstrapInv(N = 5, filterSet = 'sub')
# Stats are output to a dataset, with variables per iteration
data.stats
# Similarly to plotSpectra(), there is a plotUncertainties() function for this data
data.plotUncertainties()
# As before, an object can be returned for further plotting options
hmap = data.plotUncertainties(returnMap=True)
hmap.overlay('BLM') #.cols(1)
hmap.overlay('run').layout('BLM').cols(1)
import scooby
scooby.Report(additional=['holoviews', 'xarray'])
# tmo-dev class version
data.__version__